Fuera de la barrera,
$$ -\frac{\hbar^2}{2m}\frac{\partial ^2 \psi(x)}{\partial x^2} -E \psi(x)=0 $$de manera que las soluciones serán de la forma $A e^{{\bf i}k x}+ B e^{-{\bf i}k x}$, correspondiendo a una energía $E=\frac{\hbar^2 k^2}{2m}$. Si estamos interesados en un flujo de partículas que vienen desde $x\ll 0$, esperamos que la forma asintótica de la solución en la región $x\gg a$ sea $$ \psi(x)=A_{r} e^{{\bf i}k (x-a)} $$
En la zona de la barrera,
$$ -\frac{\hbar^2}{2m}\frac{\partial ^2 \psi(x)}{\partial x^2} + (V_0 \frac{x}{a}-E) \psi(x)=0 $$con la condición de que $\psi(a)=A_{r}$ y $\psi^{'}(a)={\bf i} k A_{r}$.
Para resolver la ecuación diferencia, conviene hacer el cambio de variables $x \rightarrow z= (x-x_0)/\alpha$, con
$$x_0=\frac{a k^2}{k_0^2}$$$$\alpha=\sqrt[3]{a/k_0^2}$$$$k_0= \sqrt{\frac{2 m V_0}{\hbar^2}}$$de manera que $\psi(x)=\varphi((x-x_0)/\alpha)=\varphi(z)$.
Aquí, $x_0$ corresponde al punto de retorno clásico de la partícula en el potencial para la energía dada, suponiendo que $a\rightarrow \infty$. Si $0<x_0<a$ diremos que el punto de retorno es real, mientras que si está fuera de la barrera tendremos un punto de retorno clásico virtual. Por otro lado, $\alpha$ una escala de longitud efectiva, que depende del potencial. Observamos que $(a/\alpha)^3=\frac{2 m a^2 V_0}{\hbar^2}=a^2 k_0^2$, de manera que esta cantidad será pequeña si $V_0$ es pequeño frente a la energía cinética mínima de un estado concentrado en la barrera.
Con este cambio, la ecuación diferencial se reduce a $$ \frac{\partial ^2 \varphi(z)}{\partial z^2} - z \varphi(z)=0 $$ que se conoce como Ecuación de Airy. Las soluciones son combinaciones lineales de las Funciones de Airy $Ai(z)$ y $Bi(z)$. Estas funciones pueden escribirse también como $$ Ai(z)=M(z)\cos(\theta(z)) $$
$$ Bi(z)=M(z)\sin(\theta(z)) $$con $M(z)$ y $\theta(z)$ dos funciones monótonas. Con la normalización usual, estas funciones tienen los siguientes límites:
Función | $z\rightarrow 0$ | $z\rightarrow \infty$ | $z\rightarrow -\infty$ |
---|---|---|---|
$M(z)$ | $\frac{2}{3^{2/3} \Gamma \left(\frac{2}{3}\right)}(1+\frac{(\frac{3}{2})^{1/3}\sqrt{\pi}}{\Gamma(1/6)}z)$ | $\frac{1}{\sqrt{\pi}z^{1/4}}e^{\frac{2}{3}z^{3/2}}$ | $\frac{1}{\sqrt{\pi}(-z)^{1/4}}$ |
$\theta(z)$ | $\frac{\pi}{6}-\frac{ 3^{5/6}\Gamma(2/3)}{2 \Gamma(1/3)}z$ | $\frac{1}{2}e^{-4/3 z^{3/2}}$ | $\frac{\pi}{4}+\frac{2}{3}\sqrt{-z^3}$ |
Derivando la identidad $\tan(\theta(z))=\frac{Ai(z)}{Bi(z)}$ respecto de $z$, encontramos además que $$ \theta^{'}(z)=\frac{-W(z) \cos^2{\theta}}{Bi(z)^2}=\frac{-W(z)}{M(z)^2}\,, $$ donde $W(z)=Ai(z)Bi^{'}(z)-Ai^{'}(z)Bi(z)=\det(^{Ai(z)}_{Ai^{'}(z)}\;^{Bi(z)}_{Bi^{'}(z)})$ es el Wronskiano asociado a las soluciones. De la antisimetría, la bilinealidad del determinante, y de la ecuación diferencial que define a $Ai$ y $Bi$, $\frac{dW}{dz}=0$ y por lo tanto $W(z)=W(0)=\frac{1}{\pi}$. De esta manera, $$ \theta^{'}(z)=\frac{-1}{\pi M^2(z)} $$
y
$$ Ai^{'}(z)= M(z)\left(\frac{M^{'}(z)}{M(z)}\sin(\theta(z)) - \frac{\cos(\theta(z))}{\pi M^2(z)}\right) $$$$ Bi^{'}(z)= M(z)\left(\frac{M^{'}(z)}{M(z)}\cos(\theta(z)) + \frac{\sin(\theta(z))}{\pi M^2(z)}\right) $$
In [105]:
%matplotlib inline
import numpy as np
import scipy.special as scpsp
import matplotlib.pyplot as plt
xvals = np.linspace(-20,1,500)
ai,aip,bi,bip=scpsp.airy(xvals)
m=np.sqrt(ai**2+bi**2)
plt.plot(xvals,ai,label="Ai(x)")
plt.plot(xvals,bi,label="Bi(x)")
plt.plot(xvals,m,label="M(x)")
plt.legend()
plt.show()
xvals = np.linspace(-2,2,500)
ai,aip,bi,bip = scpsp.airy(xvals)
m = np.sqrt(ai**2+bi**2)
theta = np.arctan(ai/bi)
theta = np.array([t+ 3.1415926 if t<0 else t for t in theta])
thetaapp0 = np.array([ np.pi/6.-.63134 * t for t in xvals])
thetaappp = np.array([ .5*np.exp(-4/3.*t**(3/2)) for t in xvals])
thetaappm = np.array([ np.pi/4 + 2./3. * np.sqrt(-t**3) for t in xvals])
plt.title("$\\theta(x)$")
plt.plot(xvals,theta,label="$exacta$")
plt.plot(xvals,thetaapp0,label="$x \\approx 0$" )
plt.scatter(xvals,thetaappp,color="g",marker="+",label="$x\\rightarrow \\infty$" )
plt.scatter(xvals,thetaappm,color="g",marker="*",label="$x\\rightarrow -\\infty$" )
plt.legend()
plt.show()
m0 = np.array([.71 +.1957*t for t in xvals])
mp = np.array([ np.exp(2/3 *t**(3/2))/t**(1./4.)/np.sqrt(np.pi) for t in xvals])
mm = np.array([1./(-t)**(1./4.)/np.sqrt(np.pi) for t in xvals])
plt.title("M(x)")
plt.plot(xvals,m,label="exacta")
plt.plot(xvals,m0,label="$x \\approx 0$" )
plt.scatter(xvals,mp,color="g",marker="+",label="$x\\rightarrow \\infty$" )
plt.scatter(xvals,mm,color="g",marker="*",label="$x\\rightarrow -\\infty$" )
plt.legend(loc=0)
plt.show()
Conviene entonces escribir la solución de la ecuación de onda con energía $E=\frac{\hbar^2 k^2}{2 m}$ como
$$ \psi(x)=\left\{ \begin{array}{l r} A_l e^{{\bf i}k x}+B_l e^{-{\bf i}k x} & x <a \\ \left.u Ai(z)+v Bi(z)\right|_{z=\frac{x}{\alpha}-(k \alpha)^2} & 0 \leq x\leq a\\ A_r e^{{\bf i}k (x-a)} + B_r e^{-{\bf i}k (x-a)} & x>a \end{array} \right. $$Escribimos las condiciones de continuidad en $x=0$
$$ \left(^{Ai(z_l)}_{Ai(z_l)/\alpha}\;^{Bi(z_l)}_{Bi(z_l)/\alpha} \right)\left(^{u}_{v} \right)= \left(^{1}_{{\bf i k}}\;^{1}_{-{\bf i k}} \right) \left(^{A_l}_{B_l} \right) $$y en $x=a$,
$$ \left(^{Ai(z_r)}_{Ai(z_r)/\alpha}\;^{Bi(z_r)}_{Bi(z_r)/\alpha} \right)\left(^{u}_{v} \right)= \left(^{1}_{{\bf i k}}\;^{1}_{-{\bf i k}} \right) \left(^{A_r}_{B_r} \right) $$donde $z_l=\frac{a}{\alpha}$ y $z_r=\frac{a}{\alpha}-(\alpha k)^2$. De esta manera,
\begin{eqnarray*} \left(^{A_l}_{B_l} \right) &=&\left(^{1}_{{\bf i} k}\;^{1}_{-{\bf i }k} \right)^{-1}\left(^{Ai(z_l)}_{Ai^{'}(z_l)/\alpha}\;^{Bi(z_l)}_{Bi^{'}(z_l)/\alpha} \right) \left(^{Ai(z_r)}_{Ai^{'}(z_r)/\alpha}\;^{Bi(z_r)}_{Bi^{'}(z_r)/\alpha} \right)^{-1} \left(^{1}_{{\bf i }k}\;^{1}_{-{\bf i }k} \right) \left(^{A_r}_{B_r} \right)\\ &=&\frac{\alpha^2}{2{\bf i}k \alpha W(z_r)}\left(^{{\bf i} k}_{{\bf i }k}\;^{1}_{-1} \right)\left(^{Ai(z_l)}_{Ai^{'}(z_l)/\alpha}\;^{Bi(z_l)}_{Bi^{'}(z_l)/\alpha} \right) \left(^{Bi^{'}(z_r)/\alpha}_{-Ai^{'}(z_r)/\alpha} \; ^{-Bi(z_r)}_{Ai(z_r)} \right) \left(^{1}_{{\bf i }k}\;^{1}_{-{\bf i }k} \right) \left(^{A_r}_{B_r} \right)\\ &=&\frac{\pi}{2{\bf i}k \alpha} \left(^{{\bf i} k \alpha}_{{\bf i }k \alpha}\;^{1}_{-1} \right) \left(^{Ai(z_l)}_{Ai^{'}(z_l)}\;^{Bi(z_l)}_{Bi^{'}(z_l)} \right) \left(^{Bi^{'}(z_r)}_{-Ai^{'}(z_r)} \; ^{-Bi(z_r)}_{Ai(z_r)} \right) \left(^{1}_{{\bf i }k \alpha}\;^{1}_{-{\bf i }k \alpha} \right) \left(^{A_r}_{B_r} \right) \end{eqnarray*}con $W(z)=\det\left(^{Ai(z)}_{Ai^{'}(z)}\;^{Bi(z)}_{Bi^{'}(z)} \right)=W(0)=\frac{1}{\pi}$ es el Wronskiano asociado a las funciones de Airy, y donde usamos que $$ \left(^{a}_{c}\;^{b}_d\right)^{-1}=\frac{1}{\det \left(^{a}_{c}\;^{b}_d\right)}\left(^{d}_{-c}\;^{-b}_a\right) $$
Consideremos ahora el coeficiente de transmisión $T$ para una partícula que se mueve desde la izquierda. En tal caso, $B_r=0$, y $T=\frac{|A_r|^2}{|A_l|^2}$.
Consideremos varios límites.
Si $k \alpha\ll 1$, \begin{eqnarray*} \left(^{A_l}_{B_l} \right) &\approx& \frac{\pi}{2{\bf i}k \alpha} \left(^{0}_{0}\;^{1}_{-1} \right) \left(^{Ai(0)}_{Ai^{'}(0)}\;^{Bi(0)}_{Bi^{'}(0)} \right) \left(^{Bi^{'}(a/\alpha)}_{-Ai^{'}(a/\alpha)} \; ^{-Bi(a/\alpha)}_{Ai(a/\alpha)} \right) \left(^{1}_{0}\;^{1}_{0} \right) \left(^{1}_{0} \right)\\ &=& \frac{\pi Bi^{'}(a/\alpha)}{2{\bf i}k \alpha} \left(^{Ai^{'}(0)}_{-Ai^{'}(0)}\;^{Bi^{'}(0)} _{-Bi^{'}(0)} \right) \left(^{1}_{-Ai^{'}(a/\alpha)/Bi^{'}(a/\alpha)} \right) \\ &=& (Ai^{'}(0)Bi^{'}(a/\alpha)-Bi^{'}(0) Ai^{'}(a/\alpha)) \frac{\pi }{2{\bf i}k \alpha} \left(^{1}_{-1}\right) \end{eqnarray*}
Si $a/\alpha \ll 1$ (barrera delgada), $Ai^{'}(a/\alpha)\rightarrow \frac{1}{3^{1/3}\Gamma(1/3)}\approx -0.26$, $Bi^{'}(a/\alpha)\rightarrow \frac{3^{1/6}}{\Gamma(1/3)}\approx 0.45\, ,$ y $$ \left(^{A_l}_{B_l} \right) \approx {\bf i}\frac{a^2}{4 k \alpha^3}\left(^{1}_{-1} \right)={\bf i}\frac{a^2 k_0^2}{4 k a}\left(^{1}_{-1} \right) $$
de donde obtenemos, para el coeficiente de transmisión, $T\approx \frac{16 k^2 /k_0^2 }{a^2 k_0^2}$
Por otro lado, si la barrera es ancha, $Ai^{'}(z_r) \rightarrow 0$, $Bi^{'}(z_r)\rightarrow \frac{e^{-2/3 (a/\alpha-(k \alpha)^2)^{3/2} }}{\sqrt{\pi}\sqrt[4]{a/\alpha}}$ y
\begin{eqnarray*} \left(^{A_l}_{B_l} \right) &\approx& {\bf i} (0.23 \frac{e^{2/3 (a/\alpha-(k \alpha)^2)^{3/2} }}{\sqrt[4]{a/\alpha}k \alpha}) \left(^{1}_{-1}\right) \end{eqnarray*}de manera que $T\approx 18.9 (k \alpha)^2 \sqrt{a/\alpha} e^{-4/3 (a/\alpha-(k \alpha)^2)^{3/2}}$
Para altas energías, $z_l < z_r \ll 1$, por lo que podemos usar la aproximación asintótica para las funciones de Airy:
\begin{eqnarray*} Ai(z)&\rightarrow &\frac{1}{\sqrt{\pi}(-z)^{-1/4}}\sin(\frac{\pi}{4}+\frac{2}{3}\sqrt{-z^3})\\ Bi(z)&\rightarrow &\frac{1}{\sqrt{\pi}(-z)^{-1/4}}\cos(\frac{\pi}{4}+\frac{2}{3}\sqrt{-z^3})\\ Ai^{'}(z)&\rightarrow &- (-z)^{1/2} Bi(z) \\ Bi^{'}(z)&\rightarrow & (-z)^{1/2} Ai(z) \end{eqnarray*}De esta manera, definiedo $\theta_{^l_r}=\frac{\pi}{4}+\frac{2}{3}\sqrt{-z_{^l_r}}$, podemos escribir la ecuación matricial como
\begin{eqnarray*} (^{A_l}_{B_l})&\approx&\frac{1}{2{\bf i}k \alpha (z_l z_r)^{1/4}} \left(^{{\bf i} k \alpha}_{{\bf i }k\alpha}\;^{1}_{-1} \right) \left(^{1}_{0}\;^{0}_{\sqrt{-z_l}} \right) \left(^{\sin(\theta_l)}_{-\cos(\theta_l)}\;^{\cos(\theta_l)}_{\sin(\theta_l)} \right) \left(^{\sin(\theta_r)}_{\cos(\theta_r)} \; ^{-\cos(\theta_r)}_{\sin(\theta_r)} \right) \left(^{\sqrt{-z_r}}_{0}\;^{0}_{1} \right) \left(^{1}_{{\bf i }k \alpha}\;^{1}_{-{\bf i }k \alpha} \right) \left(^{A_r}_{B_r} \right)\\ &\approx& \frac{1}{2{\bf i}k \alpha (z_l z_r)^{1/4}} \left(^{-{\bf i} k \alpha}_{-{\bf i }k\alpha}\;^{-\sqrt{-z_l}}_{\sqrt{-z_l}} \right) \left(^{\cos(\theta_r-\theta_l)}_{-\sin(\theta_r-\theta_l)}\;^{\sin(\theta_r-\theta_l)}_{\cos(\theta_r-\theta_l)} \right) \left(^{\sqrt{-z_r}}_{{\bf i }k \alpha}\;^{\sqrt{-z_r}}_{-{\bf i }k \alpha} \right) \left(^{A_r}_{B_r} \right)\\ &\approx& \frac{\sqrt{k \alpha}}{2{\bf i} \sqrt[4]{1 - k_0^2/k^2}} \left(^{-{\bf i} }_{-{\bf i }}\;^{-1}_{1} \right) \left(^{\cos(\theta_r-\theta_l)}_{-\sin(\theta_r-\theta_l)}\;^{\sin(\theta_r-\theta_l)}_{\cos(\theta_r-\theta_l)} \right) \left(^{\sqrt{1 - k_0^2/k^2}}_{{\bf i }}\;^{\sqrt{1 - k_0^2/k^2}}_{-{\bf i }} \right) \left(^{A_r}_{B_r} \right)\\ &\approx& \frac{-\sqrt{k \alpha}}{2 \sqrt[4]{1 - k_0^2/k^2}} \left(^{e^{{\bf i }(\theta_r-\theta_l)}(\sqrt{1 - k_0^2/k^2}+1)}_{e^{{\bf i }(\theta_r-\theta_l)}(\sqrt{1 - k_0^2/k^2}-1)} \;\;^{e^{{\bf i }(\theta_r-\theta_l)}(\sqrt{1 - k_0^2/k^2}-1)}_{e^{{\bf i }(\theta_r-\theta_l)}(\sqrt{1 - k_0^2/k^2}+1)}\right) \left(^{A_r}_{B_r} \right) \end{eqnarray*}Finalmente, $|A_l|^2 \approx |A_r|^2\frac{k \alpha}{4 \sqrt{1-k_0^2/k^2}}(\sqrt{1 + k_0^2/k^2}-1)^2$, $|B_l|^2 \approx |A_r|^2\frac{k \alpha}{4 \sqrt{1-k_0^2/k^2}}(\sqrt{1 - k_0^2/k^2}-1)^2$ y por lo tanto, $$ R= \frac{|B_l|^2}{|A_l|^2}\approx \frac{(\sqrt{1 - k_0^2/k^2}-1)^2}{(\sqrt{1 - k_0^2/k^2}+1)^2} \approx \frac{k_0^4}{16 k^4} $$
In [1]:
%matplotlib inline
import numpy as np
import scipy.special as scpsp
import matplotlib.pyplot as plt
def theta(z):
ail = scpsp.airy(z)
ml = np.sqrt(ail[0]**2 + ail[2]**2)
thl = np.arctan2(ail[0],ail[2])
mllnp = (ail[0]*ail[1]+ail[2]*ail[3])/ml**2
thlp = 1./(np.pi*ml**2)
return thl
def coeftransm(k=0.,a = 10.**(1./3.),alpha=1.):
zl = -(k * alpha)**2
zr = a/alpha + zl
ail = scpsp.airy(zl)
ml = np.sqrt(ail[0]**2 + ail[2]**2)
thl = np.arctan2(ail[0],ail[2])
mllnp = (ail[0]*ail[1]+ail[2]*ail[3])/ml**2
thlp = -1./(np.pi*ml**2)
air = scpsp.airy(zr)
mr = np.sqrt(air[0]**2 + air[2]**2)
thr = np.arctan2(air[0],air[2])
mrlnp = (air[0]*air[1]+air[2]*air[3])/mr**2
thrp = -1./(np.pi*mr**2)
u = (1.j*k*alpha-mrlnp)*np.cos(thr)/thrp + np.sin(thr)
v = -(1.j*k*alpha-mrlnp)*np.sin(thr)/thrp + np.cos(thr)
al = (u*np.sin(thl)+v*np.cos(thl))*(1.j*k*alpha+mllnp)
al = al + thlp * ( u * np.cos(thl) - v * np.sin(thl) )
al = al/(2.*1.j*k*alpha)
t = mr**2/ml**2/(al*al.conj()).real
return (t,u,v)
def coeftransm_app(k=0.,a = 10.**(1./3.),alpha=1.):
t=.5
k0 = np.sqrt(a/alpha**3)
if k > k0:
t=1- (.5*k0/k)**4
elif k<k0:
if a/alpha>1:
t = 18.9
t = t * (k*alpha)**2*np.sqrt(a/alpha)
t = t * np.exp(-4./3.*(a/alpha)**1.5 * (1-(k/k0)**2)**1.5)
else:
t= 16.*(k/k0)**2/(a*k0)**2
if t>1:
t=np.nan
return (t,1-t)
ks = np.linspace(0,4.,100)
ens = ks**2
a=1.
k0=np.sqrt(a)
ts = np.array([coeftransm(k,a=a)[0].real for k in ks])
plt.plot(ks,ts,label="$a=\\alpha$")
tsapp = np.array([coeftransm_app(k,a=a)[0] for k in ks])
plt.scatter(ks,tsapp,label="$a= \\alpha$ (approx)",c="b")
plt.rc('text', usetex=False)
plt.xlabel("k/k0")
plt.ylabel("T")
plt.title("Barrera intermedia")
plt.legend(loc=0)
plt.show()
ks = np.linspace(0,4.,100)
a=(10.**(1/3.))
k0=np.sqrt(a)
ts = np.array([coeftransm(k,a=a)[0].real for k in ks])
plt.plot(ks,ts,label="$a=\\sqrt[3]{10} \\alpha$")
tsapp = np.array([coeftransm_app(k,a=a)[0] for k in ks])
plt.scatter(ks,tsapp,label="$a=\\sqrt[3]{10} \\alpha$ (approx)",c="b")
plt.rc('text', usetex=False)
plt.xlabel("k/k0")
plt.ylabel("T")
plt.title("Barrera ancha")
plt.legend(loc=0)
plt.show()
ks = np.linspace(0,3.,100)
a=(.1**(1/3.))
k0=np.sqrt(a)
ts = np.array([coeftransm(k,a=a)[0].real for k in ks])
plt.plot(ks,ts,label="$a=\\sqrt[3]{.1} \\alpha$")
tsapp = np.array([coeftransm_app(k,a=a)[0] for k in ks])
plt.scatter(ks,tsapp,label="$a=\\sqrt[3]{.1} \\alpha$ (approx)",c="b")
plt.rc('text', usetex=False)
plt.title("Barrera delgada")
plt.xlabel("k/k0")
plt.ylabel("T")
plt.legend(loc=4)
plt.show()
In [45]:
%matplotlib inline
import matplotlib.pyplot as plt
plt.rc('text', usetex=False)
import numpy as np
x3s=np.linspace(-3.2,-1.4,100)
wf=.1*np.exp(-(x3s+2.4)**2/.25)*np.cos(20.*(x3s+1.8))+1.
plt.plot(x3s,wf,color="b")
x3s=np.linspace(-1.2,.7,100)
ks = np.array([ 20.*(1-x) if 0<x<1 else 20. for x in x3s])
wf=.095*np.exp(-(x3s+0.2)**2/.25)*np.cos(ks*(x3s))+1.
plt.plot(x3s,wf,color="orange")
x3s=np.linspace(1.1,2.9,100)
wf=.03*np.exp(-(x3s-1.8)**2/.25)*np.cos(20*(x3s+1.8))+1.
plt.plot(x3s,wf,color="green")
xs=np.linspace(-3,3,100)
vs=np.array([1.3*x if 0.<x<1. else 0 for x in xs])
plt.xlabel("$x/a$")
plt.ylabel("$V(x)$")
plt.text(1.2,1.2,"$V_0$")
plt.plot(xs,vs,color="red")
plt.arrow(-2.2,1.15,.5,0.,color="b",shape="right",head_width=.05,head_length=.15)
plt.arrow(.1,1.15,-.5,0.,color="orange",shape="left",head_width=.05,head_length=.15)
plt.arrow(1.7,1.15,.5,0.,color="green",shape="right",head_width=.05,head_length=.15)
plt.savefig("grafico-barrera.png")
In [ ]: